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Abstract. In this study we show how errors due to finite box size af- 
fect formation and destruction rate for haloes in cosmological N-body 
simulations. In an earlier study we gave an analytic prescription of 
finding the corrections in the mass function. Following the same ap- 
proach, in this paper we give analytical expressions for corrections in 
the formation rate, destruction rate and the rate of change in comoving 
number density, and compute their expected values for the power law 
(n = -2) and LCDM models. 
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1. Introduction 

In prevalent models of structure formation in the universe, galaxies, clusters of 
galaxies and other large scale structures are believed to have formed due to grav- 
itational amplification of small perturbations (Peebles 1980; Padmanabhan 1993; 
Peacock 1999; Padmanabhan 2002; Bernardeau et al. 2002). In cold dark matter 
dominated models these perturbations collapse hierarchically i.e., smaller struc- 
tures form first and then they merge and form larger structures. 

In the models of structure formation the main goal is to understand the growth 
of perturbations at various scales. It has been found that the linear perturbation 
theory closely follows the actual growth of perturbation at any scale as long as the 
amplitude of perturbation at that scale is small. However, the linear approximation 
breaks down once the amplitude becomes large because perturbations at various 
scales couple to each other and the system becomes nonlinear. Many approxima- 
tions have been proposed (Zel'dovich 1970; Gurbatov et al. 1989; Saichev Sz 
Shandarin 1989; Matarrese et al. 1992; Brainerd et al. 1993; Hui & Bertschinger 
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1996; Bagla k Padmanabhan 1994; Sahni k Coles 1995) to explain nonlinear grav- 
itational clustering but their validity has been limited to special cases. Cosmologi- 
cal N-body simulations are the main tools to understand the gravitational clustering 
in nonlinear regime (Efstathiou G., et al. 1985; Bertschinger 1998; Bagla k Pad- 
manabhan 1997a; Bagla 2005). 

In N-body simulations we simulate a representative region of the universe which 
is in general large but finite. This restricts us from incorporating perturbations at 
scales greater than the size of the simulation box. We also cannot consider pertur- 
bations at scales smaller than the grid length. This means that there is a natural 
truncation of power spectrum at large as well as small wave numbers. It has been 
shown in the earlier studies (Peebles 1974, 1985; Little, Weinberg k Park 1991; 
Bagla k Padmanabhan 1997b) that the truncation of power spectrum at small scales 
does not affect large scale clustering in any significant way. However, significant 
effects of the truncation of power spectrum at large scales on clustering at small 
scales have been found (Gelb k Bertschinger 1994a,b; Bagla k Ray 2005; Power 
k Knebe 2006; Bagla k Prasad 2006, hereafter BP06). 

Apart from N-body simulations, an analytical prescription given by Press k 
Schechter (1974) has been used to find mass function and related quantities. In this 
prescription the fraction of mass in the collapsed objects, having mass greater than 
a certain value, is identified with the fraction of volume in the initial density field 
which had density contrast, filtered over an appropriate scale, greater than a critical 
density contrast. The Press-Schechter formalism has been used in many studies 
to compute merging and mass function (Bond et al. 1991; Bower 1991; Lacey 
k Cole 1993, 1994; Sasaki 1994; Kitayama k Suto 1996; Cohn, Bagla k White 
2001). In many of these studies the main goal has been to test the results based on 
the Press-Schechter formalism against N-body simulations. Here our goal is to find 
the analytic corrections in the formation and destruction rate of haloes in N-body 
simulations due to finite box size of the simulation box. 

The effects of the box size on the physical quantities in N-body simulations 
have been studied before. Gelb k Bertschinger (1994a) showed that the rms fluc- 
tuation in the mass and pairwise velocity dispersion at a given scale are underesti- 
mated when the size of the simulation box is reduced. Bagla k Ray (2005) com- 
puted the finite box effects on the average two point correlation function, power 
spectrum and cumulative mass function in N-body simulations. They found that 
the required box size for simulating the LCDM model at high redshift is much 
larger than what people generally use. Power k Knebe (2006) studied the effects 
of box size on kinematical properties of cold dark matter haloes and concluded that 
the distribution of internal properties of haloes (like spin parameter) is not very 
sensitive to the box size. So far most of the studies used N-body simulations for 
drawing their conclusions. In one of our earlier studies (BP06) we proposed an 
analytic prescription for estimating the corrections due to box size in N-body sim- 
ulations for rms fluctuations in mass, two point correlation function and the mass 
function. 

The mass function is an important physical quantity in nonlinear gravitational 
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clustering. It has been used to compare various theoretical models with observa- 
tions (Bagla, Padmanabhan & Narlikar 1996; White, Efstathiou & Frenk 1993; 
Ostriker & Steinhardt 1995; Ma & Bertschinger 1994; White 2002). However, it 
does not contain all the information which is needed for various physical processes. 
For example, if we want to know the comoving number density of ionizing sources 
at any redshift then the mass function contains insufficient information since not all 
bound objects become the source of ionization, only those objects become ioniz- 
ing sources which form during a particular period of time (Chiu & Ostriker 2000). 
In this situation we need to know the formation and destruction rate of collapsed 
objects as a function of redshift. The rate with which quasars form in the early 
universe depend in an important way on the merger rate which is another manifes- 
tation of the formation rate (Carlberg 1990). Apart from these there are many other 
cases like understanding the epoch of cluster formation etc., for which formation 
and destruction rate of haloes are important (Bower 1991). In the present study 
we carry forward our program and give analytic expressions for the corrections 
in the formation and destruction rate and estimate their values for the power law 
(n = -2) and LCDM models. 

The plan of this paper is as follows. In section 2 we give the basic equations 
which are needed for our analysis. We give the expression for the clustering am- 
plitude (rms fluctuations in mass), mass function and the rate of change in the 
comoving number density due to formation and destruction of haloes. In section 
2.3 we show that the rate of change of comoving number density around a mass M 
can be written in terms of the rate at which the haloes of that mass are formed and 
the rate at which the haloes of that mass are destroyed due to merging. We present 
our results in section 3 and discuss their implications for the power law and LCDM 
models. In section 4 we summarize our results. 



In order to estimate the effects of box size on physical quantities in the linear 
regime, we use the mass variance a 2 (r) (Peebles 1993; Peacock 1998; Padmanab- 
han 2002) as the base quantity. 



On the basis of the correction in a 2 (r), we find the corrections in other physical 
quantities. In BP06, we showed that the variance in mass is suppressed at all scales 
when we reduce the size of the simulation box. For example, if its actual value at 
scale r is a^r) then in the initial conditions of a simulation we obtain d 2 (r, Lbox) 
where: 



2. Basic Equations 



2.1 Clustering amplitude 




(1) 
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a 2 (r, L Box ) = al{r)-a\ (r, L Box ) (2) 

Here Lbox is the size of the simulation box and af(r, Lbox) is the correction 
due to the finite box size. Note that here a\ (r, L Box ) is a positive quantity so 
clustering amplitude cr 2 (r) is always underestimated when we reduce the size of 
the simulation box. 



2.2 Mass function and number density 

In the Press-Schechter formalism, we consider the initial density field as Gaus- 
sian random and smooth it over a filter of size r (or mass M), then the fraction of 
mass in the collapsed objects F(> M), having mass greater than M, at the final 
epoch, can be identified with the fraction of volume in the initial density field which 
had smoothed density contrast greater than some critical density contrast 5 C which 
is computed on the basis of the spherical collapse model. 

•j2a(U) 

The comoving number density Nps(M, t)dM of objects which have mass in 
the range [M, M + dM] at time t is given by 



N PS (M,t)dM = ^-x d -^ldM 
1 ; M dM 

_ [2 P0 f 6 c (t) da(M) \ / \ 

"VvrMV ct2(M) dM ) eXP { 2^(M)) dM (4) 

Here we can use S c (t) = 5 C /D(t) where D(t) is the linear growth factor which 
depends on the cosmological model being considered and 5 c (t) is taken to be 1.68 
at the present epoch. 

2.3 The rate of change of number density 

We can find the rate of change in the comoving number number density per unit 
time for objects which have mass in the range [M, M + dM] from equation (4). 



Po 
M 



1 



D 2 (t) 
5 2 



dD(t)\ 



dt 



a 2 (M)D 2 (t) 



exp 



5 C da(M)\ 
<j 2 (M) dt ) 

2a 2 (M)D 2 {t) J 



dM (5) 



We can identify the first and second terms of the right hand side with the de- 
struction and formation rate respectively (Sasaki 1994; Kitayama & Suto 1996). 
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' dN PS {M,t) 
dt 



dM = - 



1 dD(t) 
D(t) dt 
dN Deat (M,t) 
dt 



a 2 (M)D 2 {t) 



dM + 



dN F 



N PS {M, t)dM 
, rm (M, t) 



dt 



dM (6) 



The formation rate (dNp orm (M,t)/dt)dM quantifies the change in the co- 
moving number density of objects around mass M, per unit time, due to the for- 
mation of objects in that mass range when objects of mass smaller than M merge 
together. 



f dN Form (M,t) 
V dt 



dM = 



1 dD(t) 
D(t) dt 

I P° (— 

ttM \D 4 (t) 



1 



a 2 (M)D 2 {t) 
dD(t) \ 
dt ) 



N PS (M, t)dM 

81 da(M)\ 
V(M) dt ) 



x exp 



g ^ 

2a 2 (M)D 2 (t) J 



dM 



(7) 



The destruction rate (dNp est (M,t)/dt)dM quantifies the rate of change of 
comoving number density of haloes in the mass range [M, M + dM] when the 
haloes in that mass range merge together and form bigger haloes. 



/ dN Dest (M,t) \ dM = 1 dD{t) 
V dt ) D{t) dt 



2po 

TTM 



N PS (M, t)dM 



1 dD(t)\ 



x cxp 



D 2 (t) dt 
2a 2 {M)D 2 {t) 



Sr. 



a 2 {M) dt 



dM 



(8) 



3. Box corrections 

On the basis of the correction due to finite box size in a 2 (r) (equation (2)) we 
can find the corrections in the number density, the rate of change in the number 
density, the rate of formation and the rate of destruction i.e., equation (4), (5), (7), 
and (8) respectively. 

3.1 Corrections in comoving number density 
On the basis of the Press-Schechter formalism the comoving number density of 
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objects in a mass range around a mass scale can be related to the rms fluctuations 
in the mass at that scale. Using the technique which we have developed to compute 
the corrections due finite box size in the rms fluctuations in mass (see equation (2)), 
we can find the corrections in the comoving number density also. 

The comoving number density which we expect in cosmological N-body sim- 
ulations is given by 



Nps(M , t)dM - ^_|_ (--^) „ P (-^ m ) *M 

= N PS ,o(M, t)dM - N PS ,i{M, t)dM (9) 

Here Nps,o{M, t) and Nps{M, t) are the theoretical (box size infinite) and the 
actual (box size finite) comoving number densities respectively. 

(10) 

and 



N PS ,i{M,t)dM 



a\ 



0"n 



x exp 



2D 2 (t) 



(Tr 



-3/2 



do\ 
dal 



N PSfi (M,t)dM 



(11) 



or 



Nps,!(M,t)dM w - 



da 2 Q 



Saj 
2aJ 



Va 2 , 



+ ^(^4 N PS ,o(M,t)dM if aj/al < 1 

(12) 

In equation (12) the coefficient of Npsfi{M,t) changes sign at the scale for 
which do ~ 5 c /D(t)V3 = 5 c (t)/y/3 so the number density of objects below this 
scale is overestimated and the number density above this is scale underestimated. 
This is in accordance with our earlier results (BP06). 

Figure (1) shows the theoretical comoving number density Npsfl and the actual 
comoving number density Nps, which we expect in a N-body simulation for the 
power law model (n = —2) at the present epoch (z = 0). In this case we consider a 
simulation box of size 128 h~ l Mpc and normalize the initial power spectrum such 
that the scale of non-linearity at z = is 8 /i _1 Mpc. The actual and approximate 
corrections in the comoving number density are also shown in the figure. This 
figures shows that in N-body simulations the number density of large mass haloes 
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Figure 1. This figure shows the actual comoving number density NpsfidM (solid line) and 
the comoving number density Nps dM (dashed line) which we expect in a cosmological 
N-body simulation at the present epoch (z = 0) for a power law (n = —2) model (see 
equation (9, 10)). Here we consider the size of the simulation box 128 h~ x Mpc and 
normalize the initial power spectrum such that the scale of non-linearity (r n i) at z = 
is 8 h~ 1 Mpc. In the figure, the exact (dot-dashed line) and approximate (dotted line) 
corrections in the comoving number density due to finite box size are also shown (see 
equation (11, 12)). 



is underestimated. However, it is overestimated for small mass haloes. This feature 
is more evident from the actual and approximate corrections (difference between 
the theoretical and actual values). In Figure (2) we show the theoretical and the 
actual comoving number densities and corrections due to finite box size for the 
LCDM model. In this case we consider a simulation box of size 10/i -1 Mpc and 
compute physical quantities at z = 6. From Figure (1) and Figure (2) it is clear 
that the overall trend remains the same. This result is in agreement with the results 
of our earlier study (BP06). 

3.2 Formation rate 

The rate of formation of haloes of mass M which we expect in N-body simu- 
lation is given by equation (7) 
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M [(M /h)] 



Figure 2. This figure shows the actual comoving number density Nps.u dM (solid line) 
and the comoving number density Nps dM (dashed line) which we expect in a cosmo- 
logical N-body simulation for the LCDM model at z = 6 (see equation (9, 10)). Here 
the size of the simulation is taken to be 10 h~ x Mpc. We have also shown the exact (dot- 
dashed line) and approximate (dotted line) corrections in the comoving number density due 
to finite box size (see equation (11, 12)) . 



'dN Form (M,t)\ dM _ [2 Po ( 1 dD{t)\( 51 da(MY 



dt J VvrMV^ 4 (t) dt J\ a 4 (M) dt J 

dN Formfi (M, t)\ dM _ f dN Formtl (M,t) \ dM (13) 



dt J V dt 

Here (dN Formfi (M, t)/dt)dM and (dN Formjl (M, t)/dt)dM are the theoreti- 
cal formation rate and the correction term respectively. 



dN Form , (M,t)\ m = Jl po / 1 dD(t) dM 



dt J \JirM\D 4 (t) dt 

_J^dMM)_\ f 51 \ 

o*{M) dt y eXP V 2ol(M)D*{t)) UJ 
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Figure 3. This figure shows the actual formation rate (dN Form,o / 'dt)dM (solid line) and 
the formation rate (dNp orm /dt) dM (dashed line) which we expect in a cosmological N- 
body simulation for the power law model (see equation (13, 14)). All the parameters for 
this figure are the same as for Figure 1 . We have also shown the exact (dot-dashed line) 
and approximate (dotted line) corrections in the formation rate due to finite box size (see 
equation (15, 16)). 



and 



'dN FormA (M,t) 
dt 



dM = - 



x exp 
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(15) 
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( dN Form}l (M,t) 
V dt 
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dM 



daf 
da\ 



5aj 
2al 



+ 



dN Formfi {M,t) 
dt 



I dM 
(16) 



Figure (3) and Figure (4) show the theoretical formation rate {dN Form $ / dt)dM 
and the actual formation rate (dN Form /dt)dM for the power law and LCDM mod- 
els respectively. The parameters for the power law and LCDM model are the 
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M [(M /h)] 



Figure 4. This figure shows the actual formation rate (dNp ormy0 /dt)dM (solid line) and 
the formation rate (dNp orm /dt) dM (dashed line) which we expect in a cosmological 
N-body simulation for the LCDM model (see equation (13, 14)). For this figure all the 
parameters are identical to that of Figure 2. The exact (dot-dashed line) and approximate 
(dotted line) corrections in the formation rate due to finite box size are also shown (see 
equation (15, 16)). 



same as in Figure (1) and Figure (2) respectively. Since the formation rate is 
directly proportional to the comoving number density, it follows the same trend 
as the comoving number density. However, in this case the scale at which the 
correction term changes sign is different from the scale at which the correction 
term for the comoving number density changes i.e., here it is the scale for which 
Co = 5 C /D(t)y/E = S c (t)/\/5. From these figures we see that the formation rate 
of massive haloes is suppressed. However, that of low mass haloes is enhanced 
in N-body simulations when the box size is reduced. The main reason behind the 
suppression of the formation of large mass haloes is the absence of fluctuations in 
the initial density field at large scales due to truncation of power. 

3.3 Destruction rate 

Following the same approach as we have applied for the formation rate, we can 
find the corrections due to finite box size for destruction rate i.e., equation (8) also. 
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2al ' 2D 2 {t)4 



dN Destfi (M,t) 
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(20) 



Figures (5) and (6) show the destruction rates for the power law and LCDM 
model respectively. The parameters for the power law and LCDM models are the 
same as in Figure (1) and Figure (2) respectively. In this case the correction term 
changes sign at the same scale at which the correction term for the comoving num- 
ber density changes. This is because Nps^i/Npsfl and NDest,i/^Dest,o are equal 
and so Nps,i and Nr)est,i nave the same scale of zero crossing. Here also we find 
that the destruction rate for the massive haloes are suppressed. However, they are 
enhanced for the low mass haloes when we reduce the size of the simulation box. 

3.4 Rate of change on number density 
The rate of change of the number density is defined as (see equation (7) 



N PS (M,t) 
dt 



dM 



-1 dD(t) 
D(t) dt 



1 
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a 2 (M)D 2 {t) 



dN Dest {M,t)\ dM+ fdN F 
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dt 



Np S (M, t)dM 
(M,t) 



dt 



dM (21) 
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Figure 5. This figure shows the actual destruction rate (dNoestfll l dt)dM (solid line) and 
the formation rate (dNrj est / dt) dM (dashed line) which we expect in a cosmological N- 
body simulation (see equation (17, 18)) for the power law model. All the parameters for 
this figure are identical to that of Figure 1 . We have also shown the exact (dot-dashed line) 
and approximate (dotted line) corrections in the destruction rates due to finite box size (see 
equation (19,20)). 



We have already given the corrections for the formation rate (dNp orm /dt)dM 
and the destruction rate (dNoest/ 'dt)dM in the last two sections. The correction 
in the rate of change of number density (dNps(M,t)/dt)dM can be written in 
terms of the correction in the formation and destruction rates. In Figure (7) and 
Figure (8) we show the rate of change of the comoving number density and the 
exact and approximate correction terms in it for the power law and LCDM model 
respectively. From equation (21) it is clear that for a(M) < S c /D(t) the rate 
of change of number density Nps is dominated by the formation rate, and for 
a(M) > 5 C /D(t) by the destruction rate. For any time t, we can find a mass scale 
M c for which a(M c ) = 5 C /D(t) i.e., the formation and the destruction rate are 
equal and so there is no net change in the comoving number density of objects at 
that scale. In hierarchical clustering models i.e., a(M) is a decreasing function of 
mass, comoving number density at large scales mainly changes due to the formation 
of massive haloes and at small scales due to destruction of smaller haloes. 

The rate of change in the number density is underestimated at large and small 
scales, however, it is overestimated at intermediate scales. This feature is clear 
from Figure (7) in which the actual and approximate error terms are positive at 
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Figure 6. This figure shows the actual destruction rate (dNn es t,o/dt)dM (solid line) and 
the formation rate {dN^est / dt) dM (dashed line) which we expect in cosmological a N- 
body simulations(see equation (17, 18) for the LCDM simulation. All the parameters for 
this figure are identical to that of Figure 2. The exact (dot-dashed line) and approximate 
(dotted line) corrections in the destruction rate due to finite box size are also plotted (see 
equation (19,20)). 



large and small scales but they are negative at intermediate scales and we have two 
zeroes crossing. 



4. Discussion 

In the hierarchical clustering models of structure formation formation and de- 
struction of haloes is a common process. In the present study we have shown that 
the formation and destruction rate of haloes due to gravitational clustering are af- 
fected significantly if the size of the simulation box is not sufficiently large. On 
the basis of the Press-Schechter formalism we have given the analytic expressions 
for the corrections in the comoving number density, formation rate, destruction rate 
and the rate of change in the number density of haloes at a given mass scale. We 
have considered the implications of our analysis for the power law (n = — 2) and 
LCDM models. Since the box corrections are more important for models which 
have significant power at large scales, so most of the models in which there is none 
or very less power at large scales are not affected by the size of the simulation box. 
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Figure 7. This figure shows the actual merger rate {dN P sfl/ dt)dM (solid line) and the 
formation rate (dNps/dt) dM (dashed line) which we expect in a cosmological N-body 
simulation for the power law (n = —2) (see equation (21)). All the parameters for this 
figure are the same as for Figure 1 . The exact (dot-dashed line) and approximate (dotted 
line) corrections in the rate of change in the comoving number density are also plotted. 



However, models in which there is a lot of power at large scales (n is large and 
negative) the box effects can be quite large. In both the cases i.e., power law and 
LCDM models, the scales at which we have shown the corrections are far below 
the size of the simulation box. 

The main conclusions of the present study are as follows: 

• If the size of the simulation box in N4^ody simulations is not large enough 
then the clustering amplitude is underestimated at all scales. 

• At any given time, there is a scale above which merging is dominated by the 
formation rate and below which it is dominated by the destruction rate. 

• The comoving number density of haloes is underestimated at large scales and 
overestimated at small scales when we reduce the size of the simulation box. 

• The formation and destruction rate also get modified by reducing the size of 
the simulation box. Particularly, they are underestimated at large scales and 
overestimated at small scales. 
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Figure 8. This figure shows the actual merger rate {dNpsfi/ dt)dM (solid line) and the 
formation rate (dNps /dt) dM (dashed line) which we expect in a cosmological N-body 
simulation for the LCDM model. All the parameters for this figure are the same as for 
Figure 2. The exact (dot-dashed line) and approximate (dotted line) corrections in the rate 
of change in the comoving number density are also plotted. 



• The suppression of the formation rate as well as destruction rate at large 
scales is mainly due to the absence of fluctuations in the initial density field 
at large scales due to limitation of box size. 

• In N-body simulations which have small power at large scales the corrections 
due to box size can be ignored. 

Observations suggest that cosmological perturbations were present at all scales 
in the initial density field that have been probed. Particularly in dark matter models 
(Diemand et al 2005; Diemand et al 2006) the index of the power spectrum at 
small scales becomes close to —3. So in order to simulate these models one has 
to be very careful in choosing the size of the simulation box. This is because if 
the index of power spectrum is close to —3 at the box scale also or if the scale 
of non-linearity is close to the size of the simulation box then the results can be 
significantly affected by the finite box effects. 

In this and Bagla & Prasad (2006) we have presented our analytic results for 
taking into account the effects of finite box size on the mass function, formation and 
destruction rate and other physical quantities. In the next paper of this series we 
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will present a detailed comparison of our analytic results with a set of cosmological 
N-Body simulations of different box size. 
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